In this section, we define a static structural model for cereals.
The paragraphs are organized so that it is easy to see how the model is built and how we can play with parameters.
Some 3D plots are interactive, don’t hesitate to rotate them.
3.1 Create a parametric leaf
A parametric leaf is here defined by sample points (here 12), described by two sets of coordinates:
- \((x,y)\): coordinates for the midrib in a vertical plane, give the curvature of the leaf;
- \((s,r)\): curviliear abcissa (i.e. abscissa along midrib) and relative leaf width along leaf shape.
Note that \(s\) can be expressed as a function of \(x\) and \(y\), at point \(p \neq 0\), as follows: $ s_p(x,y) = $.
The origin represents the leaf base.
Code
## Imports# from installed packagesimport numpy as npimport matplotlib.pyplot as pltfrom heapq import*from scipy.interpolate import splprep, splevfrom scipy.integrate import simps, trapzfrom openalea.plantgl.allimport Vector3# from ./srcfrom cereals_leaf import leaf_shape_perez, sr_prevot, parametric_leaf# or# from simple_maize import leaf_shape_perez, sr_prevot, parametric_leaf# from fitting import leaf_shape_perezfrom generator import curvilinear_abscissefrom fitting import fit2, fit3, simplifyfrom simplification import distance, cost## Code for generating a parametric leaf for a cereal pl=parametric_leaf(nb_segment=10, insertion_angle=40, scurv=0.7, curvature=70, alpha=-2.3)fig, (ax1, ax2) = plt.subplots(nrows=2) # fig.suptitle('Parametric leaf')ax1.plot(pl[0], pl[1], '.-', c="green")ax1.set_xlabel('x')ax1.set_ylabel('y')ax1.set_title("Coordinates of the midrib in a vertical plane")ax2.plot(pl[2], pl[3], '.-', c="green")ax2.plot(pl[2], -pl[3], '.-', c="green")ax2.plot(np.arange(0,1.1,0.1), np.zeros(11), c="green", ls="dashed")ax2.set_xlabel('s')ax2.set_ylabel('r')ax2.set_title("Flattened leaf shape")plt.subplots_adjust(hspace=0.5)plt.show()
Code
## Imports# from installed packagesfrom mpl_toolkits.mplot3d import Axes3Dfrom scipy.interpolate import interp2dimport matplotlib.tri as mtri# from./srcfrom fitting import leaf_to_mesh_2d## Code for representing the parametric leaf in 3D (ignore excess lines)x=pl[0]y=pl[1]s=pl[2]r=pl[3]pts,ind=leaf_to_mesh_2d(x, y, r)xs=[pt[0] for pt in pts]ys=[pt[1] for pt in pts]zs=[pt[2] for pt in pts]X,Y=np.meshgrid(xs, ys)tri=mtri.Triangulation(xs, ys)fig=plt.figure()ax=fig.add_subplot(111, projection='3d')ax.scatter(xs,ys,zs,c="green")ax.plot(xs,ys,zs,c="green")ax.plot(xs,np.zeros(len(ys)),zs,c="green",ls="dashed")ax.set_title("3D representation of a leaf shape with curvature")plt.show()
Please ignore excess straight lines.
3.2 Generate leaf azimuth series
Leaves grow around an axis with a given phyllotaxy, represented here as an angle from leaf to leaf.
Code
## Imports# from installed packages # from itertools import cycle# from ./srcfrom plant_design import leaf_azimuth## Code for generating leaf azimuths series for a given phyllotaxynb_phy=10phyllotactic_angle=137spiral=Truephyllotactic_deviation=0plant_orientation=0la=leaf_azimuth(size=nb_phy, phyllotactic_angle=phyllotactic_angle, phyllotactic_deviation=phyllotactic_deviation, plant_orientation=plant_orientation, spiral=spiral)x=np.cos(la*np.pi/180)y=np.sin(la*np.pi/180)z=np.linspace(1,len(la)+1,len(la))fig,ax=plt.subplots(subplot_kw=dict(projection='3d'))for i,a inenumerate(la): ax.plot(np.linspace(0,x[i],2), np.linspace(0,y[i],2), [z[i],z[i]], c="green")ax.plot([0,0], [0,0], [0,z[-1]], c="green")ax.set_title("3D representation of phyllotaxy")plt.show()
3.3 Manage internode lengths
Internodes on an axis have varying lengths. The repartition of their lengths along the axis can be approximated with a geometric model.
Code
def geometric_dist(height, nb_phy, q=1):""" returns distances between individual leaves along a geometric model """if q==1: u0=float(height)/nb_phyelse: u0=height*(1.-q)/(1.-q**(nb_phy+1))return [u0*q**i for i inrange(nb_phy)]
Code
## Code for applying lengths to internodes according to a geometric modelplant_height=15q=1.5x=np.cos(la*np.pi/180)y=np.sin(la*np.pi/180)z=geometric_dist(height=plant_height, nb_phy=nb_phy, q=q)fig,ax=plt.subplots(subplot_kw=dict(projection='3d'))for i,a inenumerate(la): ax.plot(np.linspace(0,x[i],2), np.linspace(0,y[i],2), [z[i],z[i]], c="green")ax.plot([0,0], [0,0], [0,z[-1]], c="green")ax.set_title("3D representation of the repartition of internode length along the stem")plt.show()
3.4 Manage leaf lengths as a function of height
Leaves that grow on an axis have varying lengths. The repartition of their lengths along the axis can be approximated with a bell shaped model.
Code
def bell_shaped_dist(max_leaf_length, nb_phy, rmax=.7, skew=0.15):""" returns leaf area of individual leaves along bell shaped model """ k =-np.log(skew) * rmax r = np.linspace(1./ nb_phy, 1, nb_phy) relative_length = np.exp(-k / rmax * (2* (r - rmax) **2+ (r - rmax) **3))# leaf_length = relative_length / relative_length.sum() * max_leaf_length leaf_length = relative_length * max_leaf_lengthreturn leaf_length.tolist()
Code
## Code for applying lengths to leaves according to a bell shaped modelmax_leaf_length=50bsd=bell_shaped_dist(max_leaf_length=max_leaf_length, nb_phy=nb_phy, rmax=.7, skew=0.15)x=np.cos(la*np.pi/180)*bsdy=np.sin(la*np.pi/180)*bsdz=geometric_dist(height=plant_height, nb_phy=nb_phy, q=q)fig, ax = plt.subplots(subplot_kw=dict(projection='3d'))for i,a inenumerate(la): ax.plot(np.linspace(0,x[i],2), np.linspace(0,y[i],2), [z[i],z[i]], c="green")ax.plot([0,0], [0,0], [0,z[-1]], c="green")ax.set_title("3D representation of the repartition of leaf length along the stem")plt.show()
3.5 Arrange a leaf to be placed along a stem with a given inclination.
A leaf must undergo a translation to be placed against the surface of the stem (tangent to it), and a rotation to be tilted like the stem (inclination = 1 if main stem).
Code
## Imports# from installed packagesfrom math import pi, cos, sin, radiansimport openalea.plantgl.allas pgl# from ./src# from cereals_leaf import arrange_leaf# or from geometry import arrange_leaf## Code for placing a leaf against a stem element (here a cylinder), with a given inclinationstem_diameter=0.5inclination=1al=arrange_leaf(leaf=pl, stem_diameter=stem_diameter, inclination=inclination, relative=True)x=al[0]y=al[1]s=al[2]r=al[3]pts,ind=leaf_to_mesh_2d(x, y, r)xs=[pt[0] for pt in pts]ys=[pt[1] for pt in pts]zs=[pt[2] for pt in pts]X,Y=np.meshgrid(xs, ys)tri=mtri.Triangulation(xs, ys)fig=plt.figure()ax=fig.add_subplot(111, projection='3d')ax.plot(xs,ys,zs,c="green")ax.plot([xs[0],xs[0]],[ys[0],-ys[0]],[0,0],c="green")ax.plot(xs,np.zeros(len(ys)),zs,c="green",ls="dashed")radius=stem_diameter/2z=np.linspace(0, zs[-1])theta=np.linspace(0, 2*np.pi)theta_grid, z_stem=np.meshgrid(theta, z)x_stem=radius*np.cos(theta_grid)y_stem=radius*np.sin(theta_grid)ax.plot_surface(x_stem, y_stem, z_stem, color="green")ax.set_title("3D representation of the placement of a leaf along a stem")plt.show()
3.6 Build the whole plant shoot in 3D, as an MTG.
An MTG is created
Code
## Imports# from installed packagesimport openalea.plantgl.allas pglfrom openalea.mtg.turtle import TurtleFramefrom openalea.mtg import MTG, fat_mtgfrom scipy.interpolate import interp1dimport pandas# from ./srcfrom geometry import slim_cylinder, stem_mesh, _is_iterable, as_tuples, addSets, leaf_mesh, compute_element, mtg_interpreter# from cereals_leaf import leaf_meshfrom geometry import CerealsTurtle, CerealsVisitorfrom fitting import leaf_element, leaf_to_mesh_2d, leaf_to_mesh, mesh4, plantgl_shape # leaf_to_mesh_new ?from plant_design import get_form_factor, blade_dimension, stem_dimensionfrom generator import majors_axes_regression, line_projection, as_leaf, as_plant, cereals as cereals_generator# not neccessary to import, but useful to see which functions are used def build_shoot(stem_radius, insertion_heights, leaf_lengths, leaf_areas, leaf_azimuths=None, leaf_shapes=None):"""create a shoot Args: stem_radius: (float) the stem radius insertion_heights: list of each leaf insertion height leaf_lengths: list of each leaf length (blade length) leaf_areas: list of each blade area collar_visible: list of each collar height or True if the collar is visible and False if it is not leaf_shapes: list of each leaf shape, if it is not known write None leaf_azimuths: list of each leaf azimuth, if it is not known write None Returns: shoot: """ ranks =range(1, len(leaf_lengths) +1) ntop =max(ranks) - np.array(ranks) +1if leaf_shapes isNone: a_leaf = parametric_leaf() leaf_shapes = [a_leaf for r in ranks]if leaf_azimuths isNone: leaf_azimuths = leaf_azimuth(len(ranks)) leaf_azimuths[1:] = np.diff(leaf_azimuths) ff = [get_form_factor(leaf) for leaf in leaf_shapes] blades = blade_dimension(area=leaf_areas, length=leaf_lengths, ntop=ntop) stem = stem_dimension(h_ins=insertion_heights, d_internode=np.array(stem_radius) *2, ntop=ntop) df = blades.merge(stem) df['leaf_azimuth'] = leaf_azimuths df['leaf_rank'] = ranks df['leaf_shape'] = [leaf_shapes[n -1] for n in df.leaf_rank]return df, cereals_generator(plant=df)
Code
def build_shoot_w_pseudo(nb_phy, plant_height, insertion_heights, leaf_lengths, leaf_areas, pseudostem_dist=1.4, stem_dist=1.2, diam_base=2.5, diam_top=1, pseudostem_height=20, leaf_azimuths=None, leaf_shapes=None, wl=0.1):"""create a shoot, with pseudostems and stems Args: stem_radius: (float) the stem radius insertion_heights: list of each leaf insertion height leaf_lengths: list of each leaf length (blade length) leaf_areas: list of each blade area collar_visible: list of each collar height or True if the collar is visible and False if it is not leaf_shapes: list of each leaf shape, if it is not known write None leaf_azimuths: list of each leaf azimuth, if it is not known write None Returns: shoot: """ ranks =range(1, len(leaf_lengths) +1) ntop =max(ranks) - np.array(ranks) +1 nb_phy =int(nb_phy)# Lejeune an Bernier formula + col = nb_young_phy =int(round((nb_phy -1.95) /1.84/1.3))# distances between leaves pseudostem = geometric_dist(pseudostem_height, nb_young_phy, pseudostem_dist) stem = geometric_dist(plant_height - pseudostem_height, nb_phy - nb_young_phy, stem_dist) internode = pseudostem + stem# stem diameters diameter = ([diam_base] * nb_young_phy + np.linspace(diam_base, diam_top, nb_phy - nb_young_phy).tolist())if leaf_shapes isNone: a_leaf = parametric_leaf() leaf_shapes = [a_leaf for r in ranks]if leaf_azimuths isNone: leaf_azimuths = leaf_azimuth(len(ranks)) leaf_azimuths[1:] = np.diff(leaf_azimuths) ff = [get_form_factor(leaf) for leaf in leaf_shapes]# blades = blade_dimension(area=leaf_areas, length=leaf_lengths, ntop=ntop) blades = blade_dimension(length=leaf_lengths, form_factor=ff, ntop=ntop, wl=wl)# stem = stem_dimension(h_ins=insertion_heights, d_internode=diameter, ntop=ntop) stem = stem_dimension(internode=internode, d_internode=diameter, ntop=ntop) df = blades.merge(stem) df['leaf_azimuth'] = leaf_azimuths df['leaf_rank'] = ranks df['leaf_shape'] = [leaf_shapes[n -1] for n in df.leaf_rank]return df, cereals_generator(plant=df)# shoot, g = build_shoot(3.0, [2,4,6,8], [2,4,6,8], [2,4,6,8])# scene, nump = build_scene(g)# display_scene(scene)
3.7 Display scenes according to different scenarii
Code
## Imports# from installed packagesfrom oawidgets.plantgl import*# Set nice color for plantsnice_green=Color3((50,100,0))
3.7.1 A single cereal
Code
## Imports# from installed packagesfrom oawidgets.plantgl import*# from ./srcfrom display import display_mtg, build_scene, display_scene# Enable plotting with PlantGL%gui qt## Code for generating a 3D cereal shoot from descritive parameters# Parameters stem_radius=1height=1500# from crop modelnb_phy=15# fixed max nb of phytomersmax_leaf_length=70insertion_angle=40scurv=0.7curvature=80phyllotactic_angle=120spiral=True# Functions callsinsertion_heights=np.array(geometric_dist(height, nb_phy, q=1.2)) # further separate stem and pseudo stem, cf simple maizeleaf_lengths=np.array(bell_shaped_dist(max_leaf_length=max_leaf_length, nb_phy=nb_phy, rmax=0.7, skew=0.15)) # plant area --> max leaf length# leaf_areas=bell_shaped_dist(plant_area=1, nb_phy=15, rmax=0.7, skew=0.15) # cf blade_dimensiona_leaf = parametric_leaf(nb_segment=10, insertion_angle=insertion_angle, scurv=scurv, curvature=curvature, alpha=-2.3)leaf_shapes = [a_leaf for l in leaf_lengths] # possible to replace leaf_length by nb_phy or...leaf_azimuths = leaf_azimuth(size=len(leaf_lengths), phyllotactic_angle=phyllotactic_angle, phyllotactic_deviation=15, plant_orientation=0, spiral=spiral)shoot, g_single = build_shoot(stem_radius=stem_radius, insertion_heights=insertion_heights, leaf_lengths=leaf_lengths, leaf_areas=None, leaf_shapes=leaf_shapes, leaf_azimuths=leaf_azimuths)# Build and display scenescene_single, nump = build_scene(g_single, leaf_material=Material(nice_green), stem_material=Material(nice_green))# display_scene(scene_single) # display in separate windowPlantGL(scene_single) # display in notebook
Code
## Imports# from installed packagesfrom oawidgets.mtg import*## Code for exploring the MTG of the generated cereal shoot# Properties on the MTG: this exclude all the topological propertiesprint(g_single.property_names())# Retrieve one property for the MTG (dict)labels = g_single.property('label')# print(labels)length = g_single.property('length')# print(length)leaf_lengths=[]leaf_ind=[]internode_lengths=[]internode_ind=[]for k,v in length.items():if k%2==0: # could have done it using labels internode_ind.append(k) internode_lengths.append(v)else: leaf_ind.append(k) leaf_lengths.append(v)width = g_single.property('shape_max_width')# print(width)leaf_widths=[]for k,v in width.items(): leaf_widths.append(v)plt.figure()plt.plot(leaf_ind, leaf_lengths, c="green", label="Leaf lengths") # == 'shape_mature_length' for final plantplt.plot(leaf_ind, leaf_widths, c="purple", label="Leaf widths")plt.plot(internode_ind, internode_lengths, c="orange", label="Internode lengths")plt.xlabel("Vertices")plt.ylabel("Lengths (cm)")plt.title("MTG properties")plt.legend()plt.show()
The leaf lengths and widths follow the bell shaped curve described before.
The first internode in the MTG actually corresponds to the pseudostem, i.e. the about 4 to 8 short first internodes that rapidly lose their leaves. The lengths following internodes follow the geometric model described before.
Code
# classes = list(set(g.class_name(vid) for vid in g.vertices() if g.class_name(vid)))# print(classes)# def vertices(g, class_name='P'):# return [vid for vid in g.vertices() if g.class_name(vid)==class_name]# vids_U = vertices(g, 'U')# plot(g_single, selection=vids_U)# plot(g_single, selection=[vid for vid in g_single.vertices() if g_single.class_name(vid)=='LeafElement'])
3.7.2 A cereal crop with variability
Code
## Imports# from installed packagesfrom random import*# from ./srcfrom stand import agronomic_plot# Enable plotting with PlantGL%gui qt## Code for generating a random cereal crop from descritive parameters with variability # Fix a seedseed(1)# Initialize the list of plantsplants_in_crop=[]# Fixed parameters for all plantslength_plot=5width_plot=5sowing_density=10plant_density=10inter_row=0.5nplants, positions, domain, domain_area, unit = agronomic_plot(length_plot, width_plot, sowing_density, plant_density, inter_row, noise=0.1)# For loop over all the plants in the cropfor n inrange(nplants):# Parameters varying among plants stem_radius=1 height=1000*(1+random()-0.5) nb_phy=15 max_leaf_length=50*(1+random()-0.5) insertion_angle=40*(1+random()-0.5) scurv=0.7*(1+random()-0.5) curvature=70*(1+random()-0.5) phyllotactic_angle=137*(1+random()-0.5) spiral=True# Functions calls insertion_heights=np.array(geometric_dist(height, nb_phy, q=1.2)) leaf_lengths=np.array(bell_shaped_dist(max_leaf_length=max_leaf_length, nb_phy=nb_phy, rmax=0.7, skew=0.15)) a_leaf = parametric_leaf(nb_segment=10, insertion_angle=insertion_angle, scurv=scurv, curvature=curvature, alpha=-2.3) leaf_shapes = [a_leaf for l in leaf_lengths] leaf_azimuths = leaf_azimuth(size=len(leaf_lengths), phyllotactic_angle=phyllotactic_angle, phyllotactic_deviation=15, plant_orientation=random()*360, spiral=spiral) shoot, g_var = build_shoot(stem_radius=stem_radius, insertion_heights=insertion_heights, leaf_lengths=leaf_lengths, leaf_areas=None, leaf_shapes=leaf_shapes, leaf_azimuths=leaf_azimuths)# Fill the list of plants plants_in_crop.append(g_var)# Build and display scenescene_var, nump = build_scene(plants_in_crop, positions, leaf_material=Material(nice_green), stem_material=Material(nice_green))# display_scene(scene_var)PlantGL(scene_var)
3.7.3 A seemingly growing plant
Code
## Inputs# from installed packagesfrom openalea.plantgl.allimport*# Enable plotting with PlantGL%gui qt## Code for generating a "growing" cereal shoot from descritive parameters# Fix a seedseed(1)# Initialize the list of plant MTGs at different stagesgrowing_plants=[]# Parameters fixing final plant or unchanged with "growth"nplants=10positions=[(x,0,0) for x inrange(-500, 500, 100)]final_height=1000final_nb_phy=2*nplantsheights=geometric_dist(final_height, final_nb_phy, 1.2)# insertion_heights=np.array(geometric_dist(height, nb_phy, q=1.2))stem_radius=2insertion_angle=30scurv=0.7curvature=100phyllotactic_angle=137spiral=True# For loop for generating plants at different stagesfor n inrange(1,nplants+1):# Parameters varying when "growing" height=10*heights[2*n-1] nb_phy=2*n max_leaf_length=5*2*n# Functions calls insertion_heights=np.array(geometric_dist(height, nb_phy, q=1.2)) leaf_lengths=np.array(bell_shaped_dist(max_leaf_length=max_leaf_length, nb_phy=nb_phy, rmax=0.7, skew=0.15)) a_leaf = parametric_leaf(nb_segment=10, insertion_angle=insertion_angle, scurv=scurv, curvature=curvature, alpha=-2.3) leaf_shapes = [a_leaf for l in leaf_lengths] leaf_azimuths = leaf_azimuth(size=len(leaf_lengths), phyllotactic_angle=phyllotactic_angle, phyllotactic_deviation=15, plant_orientation=0, spiral=spiral) shoot, g_grow = build_shoot(stem_radius=stem_radius, insertion_heights=insertion_heights, leaf_lengths=leaf_lengths, leaf_areas=None, leaf_shapes=leaf_shapes, leaf_azimuths=leaf_azimuths)# Fill the list with MTG of "growing" plant growing_plants.append(g_grow)# Build and display scenescene_grow, nump = build_scene(growing_plants[2:], positions[2:], leaf_material=Material(nice_green), stem_material=Material(nice_green))# display_scene(scene_grow)PlantGL(scene_grow)
3.7.4 An intercrop organized in rows
Code
## Inputs# from installed packagesfrom openalea.plantgl.allimport Material, Color3, Shape, Scene, Viewer, Translated, AxisRotated# Enable plotting with PlantGL%gui qt## Code for generating an intercrop from descritive parameters# Fix a seedseed(1)def plant(height, nb_phy, max_leaf_length, phyllotactic_angle, spiral):""" return the MTG of a cereal shoot generative from descriptive parameters """ stem_radius=1 insertion_angle=30 scurv=0.7 curvature=100 insertion_heights=np.array(geometric_dist(height, nb_phy, q=1.2)) leaf_lengths=np.array(bell_shaped_dist(max_leaf_length=max_leaf_length, nb_phy=nb_phy, rmax=0.7, skew=0.15)) a_leaf = parametric_leaf(nb_segment=10, insertion_angle=insertion_angle, scurv=scurv, curvature=curvature, alpha=-2.3) leaf_shapes = [a_leaf for l in leaf_lengths] leaf_azimuths = leaf_azimuth(size=len(leaf_lengths), phyllotactic_angle=phyllotactic_angle, phyllotactic_deviation=15, plant_orientation=random()*360, spiral=spiral) shoot, g = build_shoot(stem_radius=stem_radius, insertion_heights=insertion_heights, leaf_lengths=leaf_lengths, leaf_areas=None, leaf_shapes=leaf_shapes, leaf_azimuths=leaf_azimuths)return g# Organize the plant mixture in alternate rowsn_rows =10len_rows =10d_inter =70d_intra =50def plant_in_row(i):if i%(4*d_inter)==0or i%(4*d_inter)==d_inter: return plant(height=1700, nb_phy=15, max_leaf_length=100, phyllotactic_angle=137, spiral=True)else: return plant(height=900, nb_phy=20, max_leaf_length=40, phyllotactic_angle=60, spiral=True)plants_in_intercrop = [plant_in_row(x) for x inrange(0, n_rows*d_inter, d_inter) for y inrange(0, len_rows*d_intra, d_intra)]positions=[(x,y,0) for x inrange(0, n_rows*d_inter, d_inter) for y inrange(0, len_rows*d_intra, d_intra)]# Build and display scenescene_ic, nump = build_scene(plants_in_intercrop, positions, leaf_material=Material(nice_green), stem_material=Material(nice_green))# display_scene(scene_ic)PlantGL(scene_ic)
3.8 Tillering / Branching
In progress …
Code
# Copy MTG single plantg_branch = g_single.sub_mtg(g_single.root)# g_branch.display()# REARRANGE MTG WITH SCALES: PLANT, AXE, METAMER (with internode and leaf), ELEMENTS (StemElements and LeafElements)